#include "fft.h"
void FFT_WriteData(void)
{
    for (int n = 0; n < fft_N; n++)
    {
        gsta_FFT_Data[n].In = 2 +3*cos(50*2*pi*n*Ts - pi/6) + 0.8*cos(100*2*pi*n*Ts + pi/3);
        gsta_FFT_Data[n].OutRe = 0;
        gsta_FFT_Data[n].OutIm = 0;
        gsta_FFT_Data[n].OutAmp = 0;
        gsta_FFT_Data[n].OutPhase = 0;

        gsta_FFT_DataTemp[n].Re = gsta_FFT_Data[n].In;
        gsta_FFT_DataTemp[n].Im = 0;
    }
    
}

void FFT_ReadData(void)
{
    //code
    //读取输入数据
    printf("FFT输入f[n]:\n");
    for(int n=0;n<fft_N;n++)
    {
        printf("%f\t",gsta_FFT_Data[n].In);
        if((n % 8) == 7) //一行显示8个点
        {
            printf("\n");
        }       
    }
    //读取输出数据
    printf("\nFFT输出F[k]:\n");
    for(int k=0;k<fft_N;k++)
    {
        printf("F[%d]\t %fHz\t Re:%f\t Im:%f\t Amp:%f\t Pha:%f \n",k,(float)k*fs/fft_N,gsta_FFT_Data[k].OutRe,gsta_FFT_Data[k].OutIm,gsta_FFT_Data[k].OutAmp,gsta_FFT_Data[k].OutPhase);
    } 
}

//FFT
void FFT_Calculate(complex *Data,int N)
{
    if(N <= 1)
    {
        printf("数据样本不足！\n");
        return;
    }
    complex even[N/2];
    complex odd[N/2];
    for (int i = 0; i < N/2; i++)
    {
        even[i] = Data[2*i];
        odd[i] = Data[2*i+1];
    }
    //递归调用，对奇偶序列分别FFT
    FFT_Calculate(even,N/2);
    FFT_Calculate(odd,N/2);
    //对结果进行合并
    float w_N = 2*pi/N;
    complex efujkwN_Mult_oddk;
    for (int k = 0; k < N/2; k++)
    {
        efujkwN_Mult_oddk.Re = odd[k].Re*cos(w_N*k) + odd[k].Im*sin(w_N*k);
        efujkwN_Mult_oddk.Im = -odd[k].Re*sin(w_N*k) + odd[k].Im*cos(w_N*k);
        
        Data[k].Re = even[k].Re + efujkwN_Mult_oddk.Re;
        Data[k].Im = even[k].Im + efujkwN_Mult_oddk.Im;

        Data[k+N/2].Re = even[k].Re - efujkwN_Mult_oddk.Re;
        Data[k+N/2].Im = even[k].Im - efujkwN_Mult_oddk.Im;


    }
}

void FFT_Analyze(void)
{
     //直流分量
    gsta_FFT_Data[0].OutRe = gsta_FFT_DataTemp[0].Re;
    gsta_FFT_Data[0].OutIm = gsta_FFT_DataTemp[0].Im;
    gsta_FFT_Data[0].OutAmp = sqrt(gsta_FFT_Data[0].OutRe*gsta_FFT_Data[0].OutRe + gsta_FFT_Data[0].OutIm*gsta_FFT_Data[0].OutIm); // FFT_N ;       
    gsta_FFT_Data[0].OutPhase =  atan2f(gsta_FFT_Data[0].OutIm,gsta_FFT_Data[0].OutRe);
    //fs k次谐波分量
    for (int k = 1; k < fft_N; k++)
    {
        gsta_FFT_Data[k].OutRe = gsta_FFT_DataTemp[k].Re;
        gsta_FFT_Data[k].OutIm = gsta_FFT_DataTemp[k].Im;
        gsta_FFT_Data[k].OutAmp = sqrt(gsta_FFT_Data[k].OutRe*gsta_FFT_Data[k].OutRe + gsta_FFT_Data[k].OutIm*gsta_FFT_Data[k].OutIm);   
        gsta_FFT_Data[k].OutPhase =  atan2f(gsta_FFT_Data[k].OutIm,gsta_FFT_Data[k].OutRe);
    }
    

}